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We present a fully quantum mechanical treatment of the nondegenerate optical parametric os- 
cillator both below and near threshold. This is a non-equilibrium quantum system with a critical 
point phase-transition, that is also known to exhibit strong yet easily observed squeezing and quan- 
tum entanglement. Our treatment makes use of the positive P-representation and goes beyond the 
usual linearized theory. We compare our analytical results with numerical simulations and find 
excellent agreement. We also carry out a detailed comparison of our results with those obtained 
I from stochastic electrodynamics, a theory obtained by truncating the equation of motion for the 

, Wigner function, with a view to locating regions of agreement and disagreement between the two. 

' We calculate commonly used measures of quantum behavior including entanglement, squeezing and 

^ ^ ' EPR correlations as well as higher order tripartite correlations, and show how these are modified as 

^ , the critical point is approached. In general, the critical fluctuations represent an ultimate limit to 

■ the possible entanglement that can be achieved in a nondegenerate parametric oscillator. 

I \ 

. I. INTRODUCTION 

^ ' Nonlinear optical devices such as optical parametric oscillators (OPO) and optical parametric amplifiers (OPA) 

Q^ , Q have been studied in the last forty years for providing fundamental tests of quantum mechanics, as well as their 

CN ' technological applications in areas such as frequency conversion, low noise optical measurement, and cryptography. 

\ The light beams emitted by these devices are characterized by a large amount of squeezing Q, significant quantum 

■ intensity correlations 0|, and very short correlation times between the conjugate beams 01 • The entangled nature of 

■ the photons in the down converted light has been instrumental in providing experimental demonstrations^, 0,001 of 
I the original Einstein-Podolsky-Rosen paradox and other nonclassical features of quantum mechanics. In this paper, 
' we extend the usual linear theory of the nondegenerate OPO to include nonlinear effects characteristic of the onset of 

critical fluctuations near threshold, which is the physical feature that ultimately Hmits the maximum squeezing and 
entanglement available. 

As a fundamental application of these results, we point out that in 1935, Einstein, Podolsky and Rosen Q (EPR) 

■ presented their famous argument which demonstrates that local realism is inconsistent with the completeness of 
^ ■ quantum mechanics. Their argument concerned two spatially separated particles with perfectly correlated positions 
O^- and momenta, as predicted by quantum mechanics. Related correlations for quadrature phase operators have been 
^ ' studiedir n. flU fl^ and experimentally confirmed for the output fields of the nondegenerate parametric oscillator, both 

• I— I . belowUlyj and above^3| threshold. The regime of study of these correlations however has been so far confined to 
^ ' regimes of operation where the quantum fluctuations will be small so that a linearized analysis is valid. 

Closely linked with the phenomenon of EPR correlations is that of entanglement, a key feature enabling many 
potential applications in the field of quantum information. Criteria for proving entanglement using continuous variable 
(quadrature phase amplitude) measurements have been developed by Simon and Duan et al^^. Recent experiments^, 
li^lifll have measured such continuous variable entanglement but again the studies are limited to the regime of stable, 
linearizable quantum fluctuations. In this regime Gaussian statistics apply, and the criterion developed can be 
shown0| to be both a necessary and sufficient condition for entanglement in this case. 

It is known from earlier theoretical analysis [tR IT^ of the optical properties of nonlinear crystals that, in the 
linearized or Gaussian regime, a local realistic theory based on the Wigner phase space representatio n gi ves the same 
results for the correlations between signal and idler light beams produced in nonhnear crystals ^3, 0|. While this 
is also true of many correlations in second harmonic generation [l^l, there are instances where significant differences 
exist between the predictions of the two theories [lOl- Here we calculate the EPR and entanglement measures for non- 
Gaussian fields, in precisely the type of environment where non-Gaussian behavior is expected to occur experimentally 
- that is, by considering nonlinear corrections to the usual linearized approximations used to treat the OPO below 
threshold. 

In two recent papers |21J we have carried out a fully quantum mechanical analysis of nonlinear effects and critical 
fluctuations in a degenerate OPO using the positive P-representation, and have investigated the squeezing spectra 
and triple correlations in this system both analytically as well as numerically. In particular, we have shown that, in 
this case, while the full quantum theory and the semi-classical theory disagree strongly far below threshold, there 
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is a surprising agreement between the two close to the threshold where quantum fluctuations are quite intense, 
characteristic of a mixed state of Hght in this Hmit. 

The aim of the present work is to present a similar analysis for the case of a nondegenerate optical parametric 
oscillator. Both the quantum mechanical and semi-classical analyses are carried out in parallel and are compared 
with exact numerical simulations. Special attention is paid to the behavior of this system close to the critical point 
to ascertain the limits of entanglement, EPR correlations and squeezing in this regime. We find that entanglement 
is optimised just below the critical point for output mode entanglement and squeezing, while the optimum internal 
squeezing and entanglement is achieved just above threshold. 



II. HAMILTONIAN AND STOCHASTIC EQUATIONS 

We consider here the standard model for three modes coupled by a nonlinear crystal inside a Fabry-Perot inter- 
ferometer with allowance made for coherent pumping and damping due to cavity losses. This model impHes certain 
restrictions on mode-spacing for its vaHdity, since we assume only these three modes are excited. 

A. Hamiltonian 

The Heisenberg-picture Hamiltonian that describes this open system is given bv (l^l2^l2^ 

2 

i=0 

+in (^£e-"^"'al~£*e''^°'aQ^ 
2 

+ E(«^rt+atf,^ (2.1) 

i=Q 

Here £ represents the external input field at a frequency loq, with So, ai and 0,2 represent the pump, the signal 
and the idler intra-cavity modes at frequencies ujq, uji and uj2 respectively, where ujo — 101+102- The terms Ti 
represent reservoir operators and x denotes the nonlinear coupling constant due to the second order polarizability of 
the nonlinear crystal. 

This is a driven system far from thermal equilibrium, so it is not appropriate to assume a canonical ensemble. 
Instead, the density matrix must be calculated as the solution of a master equation in the Schroedinger picture. For 
simplicity, we transform to a rotating frame in which the free-field time-evolution is removed. The master equation 
for the reduced density operator, obtained after the elimination of the reservoirs using standard techniques [2J|, is 
given by 

2 

+ E 7i {2atpa\ - a\aip - pa\ai^ (2.2) 

i=0 

where 7^ are the damping rates for the mode amplitudes. For simplicity, we assume that 71 = 72 = 7 throughout this 
paper. 

To handle master equations such as this it proves convenient to transform them into c-number Fokker-Planck 
equations or equivalently into stochastic equations using operator representation theory. Here, as in our earher works, 
we use the positive P-representation for this purpose, and we also compare these results with the commonly used 
semiclassical truncation of the Wigner representation. 
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1. Classical critical point 

The classical approximation, where all fluctuations are neglected, is obtained by simply assuming that all operator 
mean values factorize. This gives us the classical nonUnear-optical equations for at = {ai) in the form: 

-Q^ = £ - 7oao - Xaia2 
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= -70;, + Xa3-.*ao (2.3) 

For small driving fields, the stable classical below-threshold solutions are — f /70 and ai = a2 = . There 
is a classical threshold or critical point aX E = Ec ~ 77o/x- Above this threshold, the driving field is clamped at 
Oio = £c/lQ^ while the signal and idler intensities increase linearly with the input field £. This paper deals with the 
near threshold and sub-threshold regime. 



B. The positive P-representation 

Using the positive P-representation[2^, we can include correlations and fiuctuations by expanding the density 
matrix describing the system in an off-diagonal coherent state basis as 

|a)((a+)*lp^(^^^+)^6^^6^+ (2.4) 



■D 
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where a = {a(),ai,a2) and a"*" = (q;q , a^", aj) are two independent triplets of complex variables. The function 
P(a,a+) can be understood as a positive phase space distribution and, by virtue of (|2.4(l . satisfies the following 
Fokker-Planck equation^3| (assuming that boundary terms vanish on partial integration): 

d 

-h— ^ [7oaJ -I- xo-t^t - £] 

+^ [71"! - Xaoc4] + [liat - xaja2] 

+ ^ [72^2 - xaoat] + [l2at - X^o ai] 

-(XQo)+ ^ +^ ^ (x«([) >-P+(a,a+,t) ■ 



9aii9a2 da{da'^ 



(2.5) 



We note here that boundary terms are found to be exponentially suppressed for large damping[2a| - ie, 7^ ^ x, 
which corresponds to typical experimental conditions for realistic OPO's in current use. 
This can equivalently be written as the following set of Ito stochastic equations [23| 

daa = [E - 7oao - xai"2) dt 

da^ = (£* — joaQ — x'^t'^t) '^^ 

dai — [—^iai+xoit<^o)dt+{xa{)f'^'^dWi 

1 /2 

dal = {-"fiol + xa2ao ) dt + (xaj) dW^ 
da2 = (-72^2 + xa^^o) dt + (xao)^^^ dW2 



1/2 

da^ = (-72"^ + XaiOio) dt + (xaj) dW^ 



(2.6) 



where 



(dWi) = {dW2) = 
{dWidW2) = (dW^dW^) = dt (2.7) 

with all other noise correlations vanishing. These equations imply that {aia]) — (rii) = when there is no driving 
field, as physically expected for a vacuum state in a normally-ordered representation. 

Numerical simulations of these stochastic trajectories confirms the assumption of asymptotically vanishing boundary 
terms for the parameters we use, as the trajectories are strongly bounded to a compact domain. At smaller damping 
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rates, it would become important to include stochastic gauge terms|20| in the equations to eliminate boundary terms, 
but this was not found to be necessary in these calculations. In other words, while boundary terms are potentially 
present, the resulting errors are expected to be of order e"'^!^ or smaller, which is completely negligible in typical 
quantum optical systems where 7 ^ X- 



C. The semi-classical theory 

We can also transcribe the master equation as a c-number phase space evolution equation using the Wigner repre- 
sentation 

Pw{oL,OL*) = ^ / d'^z xw{z.z*)e-'^ " e~'^ ° (2.8) 

where xw{z,z*), the characteristic function for the Wigner representation , is given by 

Xwiz, z*) ^ Tr (pe'"*"'+'"-") (2.9) 

This transcription is particularly useful for semi-classical treatments. 

The equation for the Wigner function for the nondegenerate parametric amplifier that corresponds to the master 
equation Ij2.2ll turns out to be 

dPw { d 

d 

d d 

d d 
+ ^ — (72^2 - X^iao) + (72^2 - xaiOo) 

2 Fi2 pa 



Q1 Q2 g: 

+70^ — 5— T + 717; — TTT + 72 



daoda^ daidal da2da2 

Q3 Q3 



X 

4 \da1da2dag ' dalda2dao 



Pw (2.10) 



If we drop the third order derivative terms, in an approximation valid in the limit of large photon number, we can 
equate the resulting truncated Fokker-Planck equation describing the evolution of the Wigner function with a set of 
Ito stochastic equations which read as follows 





= (£■ - 7oao - xaia2) dt - 


f y/^dWo 


dag 


= {£* - -^oal - xct*iCi*2) dt + y/^dW^ 


dai 


= {-jiai + xa^ao) dt + y 




da 2 


= (-7iai + X"2ao) dt + 1 




da2 


= (-7202 + X^iao) d< -1- , 


/l2dW2 


da2 


= {-120-2 + X^iag) dt + y 





The non- vanishing noise correlations are given by 

{dW,) = 

{dW^dW*) = dt- i = 0,l,2. (2.12) 

If we compare the two sets of Ito stochastic equations, namely H2.6II and II2.11|I . we notice that the main difference 
between the two is in the structure of the noise terms. While the noise terms in the positive-P equations 112. 6|l depend 
on the pumping amplitude and the nonhnear coupHng constant, those in the Wigner representation do not. In fact they 
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correspond precisely to the noise terms that one adds, in the Hnear case, in compHance with the fluctuation-dissipation 
theorem. 

In some sense, one can interpret the noise terms in the Wigner case as accounting for vacuum fluctuations. However, 
the truncated Wigner theory must be treated cautiously, since it neglects important third-order correlations which are 
not always negligible. These equations imply that {aia}) = (rii) = 1/2 when there is no driving and no coupling, as 
expected for a vacuum state in a symmetrically-ordered representation. However, a vacuum state is not obtained semi- 
classically if there is any coupling x, even with a vacuum input, which is an unphysical feature. The full Wigner theory 
has no such limitations: but it is no longer positive-deflnite, and therefore has no equivalent stochastic formulation. 



In comparing these methods, we notice that the classical equations of course give no information at all about the 
quantum fluctuations, although they give an excellent guide to the location of the critical point when the threshold 
photon number is large. 

To include quantum effects, one might imagine that a direct numerical calculation in a photon number basis would 
be useful, provided the maximum photon number was small. We note that in a three-mode system, the Hilbert space 
dimension scales as nf^^^^, while the density matrix has n^^^components provided the boson number is bounded 
by rimaxJn practice, one flnds that typical experiments have Umax — 10"^ — 10^. This implies that neither the full 
density matrix nor even the reduced wave-function in a stochastic wave-function calculation,28j can in general be 
calculated directly with current computers, for practical reasons of memory and computational time. Direct number 
state methods are also not convenient for analytical approximations. 

Other techniques involve Feynman (or related) diagrams, using a hierarchy of correlation functions [^^l- These 
methods give useful results below threshold, and have similarities to perturbation theory using stochastic methods, 
which we discuss later. The drawbacks are that these diagrammatic methods appear less systematic than phase-space 
methods, since certain classes of diagrams are discarded, and the results usually diverge at the critical point. 

Semi-classical methods overcome some of these limitations. The Wigner truncation approximation does include 
quantum fluctuation effects, but is only valid if all photon numbers are large, which is a highly questionable assumption 
in the signal/idler modes below threshold. In addition, the neglected third-order derivative terms are at best only 
suppressed by a polynomial factor, and there are no other methods to check the accuracy of the approximation. 

By comparison, the approximation of neglecting boundary terms in the positive-P equations appears well-justified 
in these calculations as long as 7 x- If necessary - that is, if unstable trajectories occur - the technique can be 
checked with the more precise stochastic gauge approach. No evidence was found that boundary terms were significant 
here, even for the relatively low damping rates used in the numerics. As one might expect, the truncated Wigner 
method gives rise to unphysical predictions at low driving field, which does not occur with the positive-P equations - 
given the parameters used here. Accordingly, we mainly focus on the positive-P phase-space method in this paper. 



In order to understand what types of calculation to carry out for this system, it is important to identify observable 
measurements that can be carried out, and relate these to operators and their correlations. 

The positive-P stochastic method directly reproduces the normally ordered correlations and moments, while the 
Wigner representation reproduces the symmetrically ordered moments. Of course, commutation relations can always 
be used to transform one type of ordering into the other. Further, we also have to distinguish between the internal 
and external operator moments, since measurements are normally performed on output fields that are external to the 
cavity. The necessary formalism for treating external field spectra was introduced and developed by Yurke_31J, and 
by Collett and GardinerJ^. 

As we shall see, there is a direct relationship between the output field spectra of a nondegenerate OPO, and 
observable criteria for EPR correlations and entanglement. 



The squeezing in terms of the intra-cavity quadrature covariances corresponds to an instantaneous measurement of 
the field moments 



D. Comparison of methods 



III. OBSERVABLE MOMENTS AND EPR SPECTRA 



A. 



Internal moments 




(3.1) 
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where we define 



^ e-'%{t) + e^'^aUt) (3.2) 



to denote internal quadrature operators. Similarly, complex quadratures jS3 are defined as: 



2 

with a normally-ordered intra-cavity variance of: 



\(n+n + ^{xr^^'-Xr^^')) , (3.3) 



= {: (t) X^" (t) :) 



4(: (xf+^/^-X^'^/^ :), (3.4) 



If such measurements were possible, they would include contributions from all frequencies. However, it is more 
typical that one has access to spectrally resolved quadrature measurements of the output fields, and these are generally 
more useful as measures of entanglement and squeezing. 

B. External spectra 

The external field measurements are obtained from the input-output relations : 



<i>7'(t) = ^27-* a,{t) - $f (t) , (3.5) 

where ^™'{t) and <i>°"*(t) are the input and output photon fields respectively, evaluated at the output-coupHng mirror, 
and aj{t) is the intra-cavity photon field. The most efficient transport of squeezing is obtained if we assume that all 
the signal losses occur through the output coupler, so that 71 = 7°"*. We will assume this to be the case for simplicity, 
though the necessary corrections[T^ for imperfect interferometers simply involve the ratio 7°"V7i- 
The measured output spectral covariance Vfj of a general quadrature 

jie out ^ e~■^e$out|^^^ ^ e'*$]°"*(t) (3.6) 

can be written as 

Vf^{iu)6{uj + uj') = (AXf °"*(w)AX^« °"*(a;')) , (3-7) 

where the ffuctuations AX^ are defined as AXj = Xj* °"* - (Xj""*), e is a phase angle related to a phase 
sensitive local oscillator measurement, and the frequency argument denotes a Fourier transform: 



/ _^e*"*Xf ™*(0 . (3.8) 
' ■/27r 



We also introduce complex quadratures and their Fourier transforms, which are useful for computational purposes: 

j^eout ^ e-*''$r*(0 + e'^*2°''*(i) 

l*°"*(w) = / _^e*"*X*™*(0 
J V 2vr 

XfSout^^^ ^ J _^g«^t^t9o«t(i) (3,9) 

The spectral quadrature operators X^ °"*(w) are not formally hermitian, except at w = 0. 
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C. Observable quadratures 

In practice, one is mostly interested in external spectral measurements taken over a long but finite interval, after a 
steady-state is achieved. For output measurements averaged over a long time T, it is the low frequency part of the 
spectrum that is the relevant quantity, as it usually determines the maximum squeezing or entanglement possible. For 
simplicity, we will focus on the w = case, where we can define observable frequency-domain quadrature operators 
as: 



Y° 



T 
2^ 



X. 



out 



(0) 



X' 



r/2 out 



(0) 



(3.10) 



Since S{u!) T/27r when mapping a Fourier transform to a Fourier series, the complex quadrature spectrum for a 
finite time interval is normalized as: 



V\0) ^(^l^X^ °"*(0)X'' °"*(0) 
In particular, the most important spectra are the unsqueezed and squeezed spectra defined by: 

y°(o) = -{ [xi + x°\ ' + [y° - 



(3.11) 



x° - x^ 



(3.12) 



In other words, the complex quadrature spectra simply correspond to simultaneous sum and difference measurements 
on the two observed output quadratures for the signal and idler, with the precise quadratures observed adjustable via 
the local oscillator phase angle 9. 

The properties of external quadratures for 7^ are experimentally important since technical noise normally 
prohibits direct quadrature measurements at w = 0. Nevertheless, even at t^; ^ the quadratures are decomposablefl^ 
into pairs of mutually commuting hermitian operators with similar properties to the intra-cavity quadrature operators, 
by using discrete sine and cosine transforms. These results therefore hold at non-zero frequencies. 

The correlations are closely related to those proposed by EPR. We will give more details in the next section, 
explaining the relationship of this type of measurement with the EPR paradox and entanglement. 



D. Stochastic mappings of operator moments 



We now wish to relate these observed operator correlations with the stochastic correlations that are used to calculate 
them via the c-number equivalences. 



1. P-representation 



In the P-representation normally-ordered operator averages directly relate to stochastic moments with respect to 
the positive P-function: 

(:X«(t)X|(t):) = (Xj(<)X|(i))p, (3.13) 
where the internal stochastic variables corresponding to the quadratures are denoted: 

X^ ^ [aje-'" + a+ e'') . (3.14) 

Also, the positive P- spectral correlations correspond to the normally ordered, time-ordered operator correlations 
of the measured fields. This leads to the following well-known result for the general squeezing spectrum, as measured 
in an external homodyne detection scheme: 

y,^(c.)<5(c. + c^0='5,,+2y9^(Al«(a;)AX^^c.'))^ . (3.15) 
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This calculation only involves the internal stochastic quadrature spectral variables, defined as: 

AX^ (cu) = J -^e-* (X^it) ~ {X^{t))p) . (3.16) 

Note that vacuum (input) field terms do not contribute directly to this spectrum, as they have a vanishing normally- 
ordered spectrum, and are not correlated with the coherent amplitudes in the positive P- representation. 



2. Wigner representation 

In the Wigner representation, on the other hand, the moments and correlations with respect to the Wigner function 
are directly related to averages of symmetrically ordered operators. It becomes necessary to rewrite the normally- 
ordered internal field averages in terms of symmetrically ordered averages using equal-time commutators. As a result 
the two spectral orderings are related by: 

(: Xf (t) X] (<):) = {Xl (t) X] (t))^ - S., . (3.17) 

Similarly, for the normally-ordered squeezing spectrum, as measured in an external homodyne detection scheme, 
one has: 

Vfjiu)6iL0 + Lo') = ^AXf ™* (uj) AX^ (cc;'))^ • (3.18) 

Here we define Fourier transforms of fiuctuations as previously, except with respect to stochastic output fields: 

X^ = e~'''$f *(t) + e'''$]°"*(t) (3.19) 

where: 

$™*(i) = V27°"*aj - $™(i) . (3.20) 

It is essential to include the vacuum field contributions from refiected input fields, as these are correlated with the 
internal Wigner amplitudes, and hence contribute significantly to the spectrum. In fact, these input fields can be 
shown to correspond directly to the noise terms in the Wigner representation stochastic equations, leading to the 
identification: 

^ = V2^f{t) , (3.21) 

where <i>™(i) is a c- number amplitude corresponding (in the Wigner representation) to the quantum vacuum input 

field, and ($f (i)$*™(t')>iv = 6(1 ~ t')/2 . 

The fundamental property of the Wigner function is that the ensemble average of any polynomial of the random 
variable a and a* weighted by the Wigner density exactly corresponds to the Hilbert-space expectation of the cor- 
responding symmetrized product of the annihilation and creation operators. Therefore, the truncated theory with a 
positive Wigner function can be viewed as equivalent to a hidden variable theory, since one can obtain quadrature 
fiuctuation predictions by following an essentially classical prescription; in which even the noise terms have a classical 
interpretation as corresponding a form of zero-point fiuctuations. This description cannot be equivalent to quantum 
mechanics in general, but may provide results which, under some circumstances, turn out to be quite similar to the 
quantum mechanical results. 



IV. EPR-CORRELATIONS AND ENTANGLEMENT 

A quantitative, experimentally testable criterion for the EPR paradox was proposed in 1989 It is important 
to understand the physical interpretation of this paradox. EPR originally assumed local realism, and claimed that 
an observation of perfectly correlated positions and momenta would imply the incompleteness of quantum mechanics. 
A modern interpretation is that one can merely deduce the inconsistency of local realism with quantum mechanical 
completeness, since local realism in Einstein's original sense is no longer widely accepted. This is a weaker paradox 
than the Bell inequality - which rules out all local realistic interpretations. However, the Bell inequality has not 
yet been violated, due to causality and/or measurement inefficiency issues (though weaker inequalities have been 
violated) . The EPR paradox with quadrature variables has the advantage that the required degree of measurement 
efficiency is readily achievable with photo-detectors, since it does not require single-photon counting. 
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A. 1989 EPR criterion: Violation of an Inferred Heisenberg Uncertainty Principle 

Consider two spatially separated subsystems at A and B. Observables Xi ("position") and Yi ("momentum") 
are defined for subsystem A, where the two operators have a commutator of Xi,Yi = 2, so that by Heisenberg's 

uncertainty principle, A^XiA^Yi > 1. Suppose that the two subsystems are partially correlated, as may occur in 
a real experiment, as opposed to the ideal correlations in the EPR ^edanfcen-experiment. One may still predict the 
result of measurement Xi, based on the result of a causally separated measurement X2 performed at B. However, 
the prediction is imperfect, and has an associated inference error. Also, for a different choice of measurement Y2 at 
B, suppose that one may predict the result of measurement Yi at A. 
We define 



Aj^fXi = J P{X2)AiXi\X2)dX2 

AlfYi - J P{Y2)A{Y,\Y2)dY2 (4.1) 



Here X2 labels all outcomes of the measurement X2 at B, and A{Xi\X2) is the standard deviation of the conditional 
distribution P{Xi\X2), where Xi is the conditional result of the measurement Xi at A, given the measurement X2 
at B. The probability P{X2) is the probability for a result X2 upon measurement of X2- 

Next, we define an inference variance A^^^Xi as the average variance of the conditional (inference) variances 

A{Xi\X2) for the prediction (inference) of the result Xi for Xi at A, conditional on a measurement X2 at B. We 
define A(Fi|l2) similarly to represent the weighted variance associated with the prediction (inference) of the result 
Yi at A, based on the result of the measurement at B. 

The 1989 inferred H.U.P. criterion to demonstrate EPR correlations in the spirit of the EPR paradox is 

A„/XiA™/Yi < 1. (4.2) 

This EPR-style criterion 114. 2|l was not given in the EPR paper, but has the useful property that it represents a 
quantitative inequality that can be experimentally satisfied, without having to construct an experimentally impossible 
state with perfect correlations, as in the original proposal. As an added advantage, the appHcation of this inequality 
to electromagnetic quadrature variables allows the use of efficient photo-detection techniques, which makes this a 
completely practical measure. 

By contrast, the violation of a Bell inequality - while having stronger consequences - is more difficult to achieve, 
owing to poor efficiencies encountered in single-particle detectors and polarizers. For either type of experiment, a 
crucial element is the causal separation of detectors. Without this, arguments using causality provide no constraints 
or inequalities at all. 



1. Linear estimate criterion 



It is not always convenient to measure each conditional distribution P{Xi\X2) and P(Yi\Y2) and its associated mean 
and variance. A simpler procedure^^ is to propose that upon a result X2 for the measurement at B the predicted 
value for the result Xi at A is given linearly by the estimate Xest = CX2 + d. The RMS error in this estimate after 
optimizing for d is 

Alf^^X, = (^2) _ ^Sof, (4.3) 

where 60 = Xi — 0X2- The best choice for c minimizes Af^j j^X and can be adjusted by experiment, or calculated as 

discussed in [ill to be c = ((^1, ^2))/A2X2 , where we define {Xi,X2) = {X1X2) - {Xi){X2). 

Generally the linear estimate will correspond not be the best estimate for the outcome at A, based on the result at 
B. Therefore generally we have Amfx^ > AinfX and Ai„/ ^lY > A„/y[i3. The observation of 

A„/,lXiA„/,lYi < 1 (4.4) 

will then also imply EPR correlations in the spirit of the EPR paradox. 
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B. An entanglement criterion based on the observation of two-mode squeezing 

Entanglement may be deduced through a whole set of criteria, of which the EPR criterion II4.2II is one[^. It is 
possible to deduce entanglement through other criteria[3 without the need to prove the strong EPR correlations. 
This has significance within quantum mechanics, but not necessarily the broader implications of the EPR criterion. 

Such entanglement criteria, derived by Duan et.al. and Simon^lJl, are based on the proof of quantum inseparability, 
where the failure of a separable density matrix 

P^Y^PrPrPR (4.5) 

R 

{J2rPr = 1) is proved. Particularly useful for our purposes is a criterion considered by Duan et al, sufficient to 
demonstrate entanglement (inseparability). We define 

SX = Xi — X2 

SY = Y1+Y2 ■ (4.6) 

Entanglement is guaranteed provided that the sum of the variances is bounded by: 

A^5X + A^5Y < 4 . (4.7) 

This observation of this entanglement criterion (I4.7|l may be identified as a "two-mode squeezing" criterion for 
entanglement, since the individual criterion 

A^6X^({X,-{Xi)-{X2-{X2})y} < 2 (4.8) 

is the criterion for the observation of a type of "two-mode squeezing". In this way we see that fields that are two-mode 
squeezed with respect to both Xi — X2 and Yi +Y2, each satisfying H4.8|l . are necessarily entangled. 



C. EPR correlations and entanglement of the parametric system 

The EPR correlations are predicted possible for the outputs of the parametric oscillator. For intracavity entangle- 
ment, we define the quadrature phase amplitudes 

Xi = {ai+a\) 

Yi = (ai - d\)/i 

X2 = {0-2 + al) 

Y2 = (a2-4)A- (4-9) 

and identify correlated observables for the oscillator, so that Xi is correlated with X2 and Yi is correlated with —Y2. 
The Heisenberg uncertainty relation for the orthogonal amplitudes of mode di is A^XiA^Yi > 1. 

As explained in the previous section, for practical reasons it is preferable to use the corresponding observable external 
quadratures defined at or near zero frequency, which are X",Y°. However, the detailed arguments only depend on 
having the commutators defined above, together with the requirement of causality - that is, the observations must 
take place with space-like separations between the two detectors over the whole observation period T. 

We calculate several types of EPR or entanglement measures. Firstly we evaluate the the 1989 inferred H. U. P. 
EPR criterion Il4.2|l but using the linear estimate form, which will allow demonstration of both entanglement and EPR 
correlations defined in the spirit of the original EPR paradox. In terms of quadrature phase amplitude measurements 
this strong EPR criterion is satisfied when 

Alf^^X^Alf^^Y" = A\X° - c.X°)A'{Y° - CyY°) < 1 (4.10) 

Now = {Xf,X^)/A'^X^ and Cy = (Yf ,Y^) / A'^Y^ will minimize^] the inference variances. Substituting for 
and Cy, we explicitly calculate 



Alf^^X" = A'X° - (X°,X°)VA2X° 



(4.11) 
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and 

Alf^^Y" = A2y° - {Y°, Y°)yA'Y° . (4.12) 

For our particular system moments such as have (ai) = (02) — ... are zero and we have symmetry between aiand 02 
modes, so that 

A^X° ^l(v + V"/'^) > 1 (4.13) 



2 



and 



(X°,X°)-i(y"-V^-/2) . (4.14) 

The linear inference EPR criterion 114. 4|l is then equivalent to: 

„ 2F°F''/2 
A^„frX'' = — ^<1. (4.15) 

This criterion is not equivalent to H4.2|l which is based on the conditionals, since the linear estimate may not be the 
best, in which case it is possible that H4.2II is satisfied while Il4.1fl|l is not, and we do not pick up EPR and entanglement 
where it exists. Nevertheless the criterion II4.15|I is sufficient to prove EPR correlations and entanglement. 

Secondly, we calculate the Duan and Simon et. al. two-mode squeezing criteria Ij4.7|l for entanglement. Written in 
terms of quadrature phase amplitude measurements this becomes 

y^/^ ^ 1 (a2(x° - X°) + A^Y" + Y°)) < 1. (4.16) 

This criterion was explicitly shown to be both sufficient and necessary for entanglement for the case of Gaussian states 
(for appropriately chosen quadratures), meaning that in this case it would pick up any entanglement present. Our 
system is not Gaussian, and while these criteria are still sufficient to imply entanglement, they may not be necessary. 

It is always the case that for ideal squeezing, both the linear EPR and the squeezed entanglement criteria are 
satisfied. Where one has additional loss, however, it is possible for the squeezed-entanglement criterion H4.16|l to 
be satisfied but not the EPR criterion ()4.inil . Such situations have been studied by Bowen et alfZl. Our situation 
could be different again, due to the fact that the underlying quantum states undergo nonhnear fiuctuations and are 
inherently non-Gaussian. 



V. BELOW-THRESHOLD INTRA-CAVITY MOMENTS 



In this section we use perturbation methods to study the nondegenerate parametric oscillator beyond the linearized 
regime both in the fully quantum mechanical approach using positive-P representation, and in the semi-classical 
approach based on the Wigner function. In the positive P case the basic quantities investigated are correlations 
involving the internal quadrature operators, mapped into stochastic variables according to 

X = (ai + a^) Y — — (ai — aj) 

X+ = (02 + a+) = 7 ("2 - at) (5.1) 

In the truncated Wigner (semi-classical) case, we have a similar set with replaced by a*. To avoid excessive 
notation we use the same symbols for the quadrature variables in the two cases, noting that in the semi-classical case, 
X+ = X* and Y+ = Y* . 

For developing a systematic perturbation procedure, it proves convenient to define 

7r = 7o/7 , fi^£/£c, ff= — ^= (5.2) 
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and to introduce the following scaled quadrature variables 



X = gX 

y = gY 



X ' 



gx+ 



y+ = gY+ . 

In terms of these new variables, and a scaled time r = jt, the equations for the quadratures become 
• Positive-P equations 



dxo 
dyo 

dx 
dy 
dx+ 
dy+ 



-Ir [xq - 2/i + {xx+ - yy^)] dr 
-Ir [yo + {xy^ + yx'^)] dr 



-a; + - {xxo + yyo) 

-y + \ i^yo - yxo) 



dT + 



dr 



-x+ + - {x+xo + y+2/o) 

-y^ + \ {x'^ya - y^xo) 



V2 

■ g 

\/2 



^xq-\- iyodwi + ^/xo- iyodw^ 
\/xo+ iyodwi - \/ xo - iyodw2 



V2 



^/xo + iyodw2 + \/xo- iyodwf 



dr — i^7= 
V2 



\/xo + iyodw2 - \/xo- iyodw^ 



where: {dwidw2) = {dwidw2) = dr . 
Semi-classical equations 



(5.3) 



(5.4) 



dxo = 
dyo 

'/,/■ = 

dy = 

dx+ = 

dy+ = 



-7r [xo -2^1+ (a;a;+ - t/j/+)] dr + \/2g^r [dwo + dw^] 
-Ir [yo + {xy~^ + ya;+)] dr - i\Pig^r \dwo - dwQ\ 



-a; + - {xxo + yyo) 

-v^\ i^yo - yxo) 



+ 2 (^^^0 + y^yo) 



-y^ + 2 {^^yo - y^xo) 



dr + g [dwi + dw2] 
dr — ig [dwi — dw2\ 

dr + g [dw2 + dwi] 
dr — ig [dw2 — dw^] . 



(5.5) 



where: {dwid,w*) = Sijdr . 

In order to solve these coupled equations systematically, we introduce a formal perturbation expansion in powers 
of 5 



En (n) 
9 xl 



n=0 



En (n) 
9 yl 



(5.6) 



n=0 



This expansion has the property that the zero-th order term corresponds to the large classical field of order l/ff in 

the unsealed quadratures, while the first order term involves the quantum fluctuation of order 1, and the higher order 
terms correspond to nonlinear corrections to the quantum fluctuations of order g and higher. 
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A. Matched power equations in Positive-P representation 



Substituting Il5.6|l in Ij5.4|l and equating like powers of g on both sides we obtain a hierarchy of stochastic equations. 
The set of equations thus obtained, if desired, can be diagrammatically analyzed using the 'stochastic diagram' 
method[23|- The zero-th order equations are 



-Ir 
-It 



(0) 



,(0) 



_2.+ (0) 



-y 



+ (0) 



2^+ {x^^'^x 



^(0)y+(0) 



(0)^+(0) 



+ (0)^ 
dT 



,,(0)^(0) 

y ■'^0 



dT 

dT 



2'" 



-(0)^(0) 



y 



dT 



dT 



dT 



(5.7) 



These equations correspond to the classical nonlinear equations for the intra-cavity quadratures expressed in terms 
of dimensionless scaled fields. Below threshold, the steady-state solution of these equations is well known and is given 
by: 



The first order equations are 



„(o) 



2n 





-JrXQ^^dT 






(1) J 






-{I- ^i)x'^^^dT + y 


/2ji.dwxi 




- {I + ^i) v'^^Ut - i 


\/%idwyi 


dx+^'^ = 




\/%LdWx2 


dy+^^^ = 


-(I + m) y+^^^dr- 


i\/%ldwy2 



(5.9) 

We have introduced new Wiener increments as dWxi(yi){T) — {dwi{T) ±(iw^(r))/V2 and dwx2(y2){'i') — {dw2{T) ± 
dw^ [t)) / ^/2, with the following correlations 



{dwxidwx2) 
{dwyidwy2) 



dT 

dT . 



(5.10) 



and all other correlations vanishing. 

The equations 115. 9|l are the ones that are normally used to predict squeezing. They are linear stochastic equations 
with non classical Gaussian white noise and, if higher-order corrections are ignored, yield an ideal squeezed state for 
the sub-harmonic quadratures together with an ideal coherent state for the pump. Further, from the structure of 
these equations, it is evident that the steady state solution for the pump field quadratures, in this order, vanish. We 
can, therefore, without loss of generality, set all odd orders of x'^\ j/q"^ for the pump, and all even orders of x\"'\ 

(n) 

Dl ' ; i = 1, 2 for the signal and idler fields respectively equal to zero. With this in mind, the second order equations 
turn out to be 



dx, 



(2) 



-It 



dT 



(5.11) 



14 



Since, in the present work, our primary interest is to calculate the first nonlinear corrections to ideal squeezed-state 
behavior, to be consistent, we need to include contributions from the third order equations as well. These equations 
are as given below 



(2) 

Vo 



dr 



1 



1 



dr 
dr 



xf^dwxi + iyl^^dwyi 



y'i^'dwxi - ix'^Q^dwyi 



1 



1 



x'^^dWa:2 + Wq'' ' dWy2 



y^^^dwx2 - iXfpdwy2 



(5.12) 



This set of equations has non-trivial noise terms as they depend on the solutions of the stochastic equations at the 
second order. 



B. Operator moments in the positive P-representation 



The set of stochastic equations equations together with the Ito rules for variable changesj23 permit computation 
of the operator moments in a straightforward manner. Apart from their intrinsic interest, they are useful in checking 
the correctness of somewhat more involved spectral calculations given later. The results obtained are summarized 
below: 



(y(l)y+(3)) 



1 -/i2 



1 +yU 



1-M 



4(1+^) (1-^^/2) 
' ^ilr ^ 7r (2-A* + /i') +4(l+/i) 

7^ + 2 (1 + At) (7r + 2(1 + n)) 

,2 



1 -Ai2 



7r 



7r 



(5.13) 



The first quantity above pertains to the depletion of the pump that supplies energy for the subharmonic mode. The 
next two quantities are the squeezed and enhanced quadratures as given by the linearized theory, while the fourth 
one is the first correction to the linearized theory. The last one is the steady state triple quadrature correlation. 
This quantity has been investigated earlier for its relevance in distinguishing quantum mechanics from a local hidden 
variable theory Q|. 

The results above yield the following expression for the steady state intra-cavity squeezed quadrature fiuctuations: 



1 



l + /i 2(1 + Ai)(l-At2) 



Ir + 



7^2 - At + M^) + 4 (1 + pi) 

(1+M) (7r + 2(l+M)) 



(5.14) 



Note that the intracavity squeezing quadrature near threshold is not perfectly squeezed, having a limiting squeez- 
ing/entanglement of 0.5, as shown in Fig (QJ. 

As might be expected, the nonlinear correction is divergent at the threshold, and needs to be handled either by 
numerical integration or a critical-point expansion. Questions relating to optimal output entanglement and squeezing 
will be treated in the next section, using frequency domain methods. 
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FIG. 1: Graph of total squeezing/entanglement moment {Y-i) vs driving field /i, using parameters of 5^ = 0.001, and 7,- = 0.5. 
This demonstrates the existence of intra-cavity squeezing and entanglement for ji <1. 



C. Matched power equations in semi-classical theory 

Using the same technique of matching the powers of g, we obtain the following set of equations in the semi-classical 
theory. The zero-th order equation are: 



dT 



j^(0) + (^(0)^+(0)+y(0)^+(0)^ 

) 



dT 
dT 



dT 



dT 
dT . 



As in the positive-P case, the steady-state solution of these equations is given by: 



2n 



,(0) 



. 



The first order equations are 



where 



dXo' = -JtXq^ dT + 2jrdWxO 

dyo^ = -Iry'o' dT + 2-^rdWyo 

dx'^^'i = - {I- ii)x^^Ut + \/2dwx-L 

dy^^^ = -{l + i^)y^^UT + V2dWyi 

= - {1 - IJ,) x+^^Ut + V2dWa^2 

dy+W = -{l+iJ,)y+<^'^UT + V2dWy2 



{dwxodwxo) = {dwyodwyo) = {dwxidwx2) = {dwyidwy2) = dT 



(5.15) 



(5.16) 



(5.17) 
(5.18) 
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with all other correlations vanishing. 

The equations above give the Hnearized theory. The first nonlinear corrections come from the next two sets of 
inations ffiven below. 



equations given below. 

The second order equations are: 



d4^) = 



dT 
dT 



dT 

dT 



The third order equations are: 



dx^P = 



-Ir 
-It 



dT 
dT 

dT 
dT 



dT 

dT 



dT 
dT 



(5.19) 



(5.20) 



D. Operator moments in semi-classical theory 

In this case, the analogues of the results in Ij5.i;ill are found to be: 



(4^^> 

iy(^)y+(^)) 
(y(l)y+(3)) 



1 

1 

1 + n 



1 



7r 



1 



7r 



2(l-;u)(l + /i) V7r + 2y 2(1 + ^)2 V7r + 2(1+At), 
7r ■^ , A* , 



4(1 - + /i)2 \-fr + 2j 2(1 - Ai)(l + yu)3 4(1 + /i)3 

1 / 7r 



7r 



7r + 2(1 + /i) 



l-M^ V7 



(5.21) 



The main difference in these calculation, compared with the positive-P results, appears in the nonlinear correction 
for the subharmonic squeezed quadrature. Up to second order in g we have 



1 



17 



10 



> 

V 




FIG. 2: Graph of second order nonlinear correction to the squeezing/entanglement moment {Yi}nl vs driving field /i, using 
parameters of = 0.001, and jr = 0.1, 1, 10. Upper lines have larger 7r values. Solid lines are the positive-P result, which 
vanish at small driving field. Dotted lines are the (less accurate) semi-classical result, which do not vanish at small driving 
field. 



1 + 2(l + /i)(l-At2) 



7^(1 + 3(1- 2(?) + 4/^(1 + /i) 



Ir ■ 



(1+M)[7r-+2(1+M)] 



(5.22) 



The similarities and disagreement between this result and the positive-P expression for the same quantity deserve 
further comments given in the concluding section. In particular, we note that, while the linear terms agree, the 
nonlinear term are not in agreement below threshold. However, just below threshold both theories give essentially 
identical nonlinear corrections. There is good agreement also in the Hmit 7,- 0. 

These comparisons are shown in Fig 

VI. SPECTRAL CORRELATIONS 

Next, we proceed to analyze spectral correlations which are of direct relevance to comparison with experiments. In 
particular, we compute the nonlinear corrections to the squeezing spectrum. 



A. Positive-P representation 



To perform calculations in the frequency domain, it proves convenient to deal directly with the Fourier transforms 



dr 



e'^^x (t) 



of the hierarchy of the stochastic equations obtained earlier. The equations thus obtained contain noise terms 



Ca:,y (^) 



with the following correlations: 



(fa m = , 



(6.1) 



In this context, for notational compactness it is useful to introduce the standard notation for convolution of two 
functions: 

[Ai.B]{n)= f ^A{n')B{n-n') . 

J V 27r 

With this in mind, the stochastic equations obtained earlier may be rewritten in the frequency domain as follows: 



• First order: 



• Second order: 



Third order: 





(0) = 




(Q) = 




m = 




(Q) = 



^xl (») 

{-iil + 1 - ^) 

{-in + i + n) 

7JI (fi) 



+ 1 - ^) 
{-in + l + n) 

74xW^f+W-yW^^+(i)] {n) 

{-in + Jr) 



{-in + Jr) 



2(-iO + 1 



2(-iO + l + ^) 



2(-ifi + l-yu) 



2 (-ifi + 1 + At) 



B. Squeezing correlation spectrum 

We now calculate the spectrum of the squeezed field, which is given by {y (Qi) y'^ (^^a))- 



{y{n,)y+{n2)) = g2^?yW(ni)y+W(02)) + 



(i^i) (^2)) + (l^i) (1^2) 



The lowest order contribution is the usual result of the Hnearized theory and given is given by 

Of + (i + m) 



In terms of the squeezing variance, this means that: 
For comparison, note that the complementary (unsqueezed) spectrum to this order is: 



4m 



+(1) 



2^1(5 (Oi +O2) 

ni + (1 - m)'" 
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Taking the next order corrections into account we find that the spectrum of the squeezed quadrature is given by 



(6.9) 



where S{n) is given by: 



s{n) = 



-2/i 



f]2 + (1+^)2 [f]2 



(02 



1 



M7r(l - M ) 

(l + ^ + 7,)(l 



(i + m)2]2 

(l-M + 7r)(l+M)-»' 

" (l-A*)[r!2 + (l-M + 7r)2] 



(l + ^)[f72 + (l + ^ + ^^)2] 

The correctness of the above expression can be checked by verifying the following equality: 



/y(l)(,)y+(3)(,)\ = 
\ / ss 



„i(ni+02)r 



2tt 



(6.10) 



(6.11) 



The corresponding external squeezing spectrum is then: 

4/i 



1 



+ (1 + /i)2 
(172 + 1 - ^2) 



[n^ + (1 



f M)2]2 
7r)(l4 



/i)-f^2 



/i7r(l - /i2) 
(l+Ai + 7r)(l 



(1-M)[f72 

M) -f^2 



(l-M + 7.)2] 



(1 + ^)[r22 + (1 + /i + 7^)2 



(6.12) 



This equation gives the complete squeezing spectrum, including all nonlinear correction to order or The 
linear part gives perfect squeezing for ^ = 1, and 17 = 0, as expected from the linear theory. The nonhnear terms give 
corrections to perfect squeezing below threshold. At zero frequency, we find that: 



l/-/2(0) 



1 



+ 



(1^ 

(i + M)^ 



M)2 

1 - 



+ 



2^2^^(2 + 7^) 



(l-A^)((l+7.)2-M2) 



(6.13) 



The resulting behavior for the optimum entanglement, which is found at zero-frequency (ignoring complications 
from technical noise), is shown in Fig J^J. We see that, as expected, the entanglement is not optimized at the critical 
point, since the nonlinear critical fiuctuations spoil this before an ideal entangled two-mode squeezed state with 
yn-/2 = Q jg achieved. Better entanglement is obtained when 7^ is reduced, as this minimizes the 'information leakage' 
in the losses of the pump mode. In this limit, the only losses are through the signal and idler output ports, which are 
needed in order to have extra-cavity measurements. 

This expression does not describe the spectrum very close to the critical point, as it diverges at the threshold. This 
region requires a different kind of scaling and is discussed later. 

The complementary or unsqueezed spectrum, for measurements of the maximum quadrature fiuctuations, is given 
by: 



y°(r!) = 1 + 



4fi 



f72- 
(^(2 



(1 

- 1 - 



2\ 



[f72 

(1- 



7r)(l- 



^) -172 



M7r(l - //2) 
(l+/i + 7,)(l 



(l+/i)[f72 + (1+^ + 7^)2 



(6.14) 



The resulting behavior for the zero-frequency critical fiuctuations is shown in Fig Near the critical point, higher 
order terms are likely to become significant. The effects of these are treated in the next section. 
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FIG. 3: Optimum squeezing with ff^ — 0.001, 7r = 10 ^, 10 ^, 10 ^, 1, 10 . Upper lines have larger values of 7r. Here, V'^^'^ < 1 
indicates squeezing and entanglement occurring at zero frequency. 



10^ 



O 2 

o— 10 
> 



10, 




0.5 



FIG. 4: Complementary (unsqueezed) spectrum with = 0.001, 7r = 10 ^,10 ^,10 ^, 1, 10 . Lower lines have larger values 

of 7r. 



We note here that in the linearized analysis, the product of these spectra corresponds to the Heisenberg uncertainty 
principal: 



1 



4m 



n■^ + {l + ^lf 



1 



4/i 



Near threshold where nonlinear effects are dominant, this relationship no longer holds. The zero-frequency nonlinear 
uncertainty product is shown in Fig ©. Just below the critical point, the nonhnear corrections apparently predict an 
uncertainty product less than unity, which is clearly the point at which the second-order perturbation method breaks 
down. An unexpected feature of these results is that for 7^ ^ 1 , the uncertainty product remains close to unity for all 
driving fields, indicating that there is a near minimum uncertainty state for low-frequency spectral measurements in 
the output fields. This does not mean that there is a minimum uncertainty state for the internal quadrature moments, 
since these are effectively integrated over all frequencies, and involve different quantum fields. 

We also investigate the behavior of the inferred Heisenberg uncertainty product, which demonstrates that there 
is an EPR paradox. In the original proposal, this uncertainty product would be zero, as the original EPR paradox 
involved perfect correlations. Instead, the minimum value of this product is determined by the nonlinear critical 
fluctuations. Due to symmetry, we only need plot the behavior of A^^^ j^X° in Fig ijSJ. 
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FIG. 5: Heisenberg uncertainty product with = 0.001, 7r ~ 10 10 , 10 ^, 1, 10 . Upper lines have larger values of 7r. 




0.5 0.6 0.7 0.8 0.9 1 



FIG. 6: Inferred quadrature uncertainty with — 0.001, 7r ~ 10 ^, 10 ^,10 ^, 1, 10 . Upper lines have larger values of 7r. 
When A^,jj^X° < 1, one may infer an EPR paradox. 

This shows quahtatively similar behavior to the entanglement measure based on squeezing, and in fact for strong 
entanglement the inferred uncertainty and squeezing measures are simply related by 

We see that near threshold, the EPR measure and squeezing entanglement measure both show the existence of a 
strongly entangled output beam, as one might expect - but the perturbation theory breaks down past the point where 
optimum entanglement is achieved, just below threshold. 

C. Triple Spectral Correlations 

Triple spectral correlations give quantum effects which distinguish very stronglyj^^ between the full quantum theory 
and the semi-classical approximation. 

Here, we calculate the internal quadrature triple spectral correlation (i (Jli) {^2) Vo (^3))- To the lowest non- 
vanishing order this is given by 

(i i^i) y+ (n^) yo (^^3)) = 9' (i^'^ i^i) (^^2) y^o^ (^3)) ■ (6.15) 
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(2) 

Substituting for yg ; have 

and using the Gaussian nature of the stochastic variables involved to factorize the fourth order correlations we 
obtain: 

x(^) in,) ^+(^) in.) yl^^ in,)) = ^,yv2.6in, + n^ + n,) _ ^^^^^^ 



V^J^/2^ S ini + 


n 


2 +f^3) 


i-in, + Jr) 


nl + (1 - m)' 




ni + ii + fi)' 



To check this result, we evaluate the steady state moment (x^^^ (r) y^^^^ (r) yg^^ (t)\ using 

\ I ss J J ^ Z-K J \/ l-K 

x(i(i)(r!i)y+(i)(r!2)yf (f^3)) • (6.19) 

and find that we obtain the same result as given earlier by direct calculations. 

This result will be compared later with the corresponding result obtained in the semi-classical theory . 

D. Comparisons with simulations 

In order to verify the accuracy of these analytic calculations, we performed extensive numerical simulations of 
the full nonlinear stochastic simulations, using a differencing technique as in earlier studies. We only calculate the 
nonlinear squeezing variance, defined as: 

Vi^) = V'l'^in) - V^^^^'I'^i^) . (6.20) 

This allows us to focus on the nonlinear corrections, which are relatively small except very near the critical threshold 
at /Lt = 1. The numerical method has the advantage that, unlike perturbation theory, it is valid at all driving fields - 
even at the critical point. 

The integration parameters used were step size c?t — 0.001, with a time- window of Tmax = 10000. The number of 
stochastic trajectories used for averaging were 2000, resulting in typical relative sampHng errors of around ±2%, as 
can be seen from the background sampling noise in some of the resulting spectra. 

Typical results are shown in Figs Ij7l8|l below, for driving fields of /i = 0.5,0.9. Note that these graphs only include 
the nonhnear corrections. Excellent agreement is found with the analytically predicted results for these values of 
driving field. 

Figure © shows results slightly closer to threshold, at /i = 0.93, which is the optimum driving field for the 
parameters chosen. 

At this point, a maximum error of around 10^** is found, due to higher order nonhnear corrections. This indicates 
that the analytic perturbation theory is able to correctly predict nonlinear effects up to the optimum squeezing point, 
but starts to diverge beyond this point. The numerical results, however, are stable throughout the critical region. To 
obtain analytic predictions in the critical region, we turn to a different asymptotic expansion in a later section. 

VII. SEMI-CLASSICAL SPECTRAL THEORY 

In this section we calculate approximate nonlinear results using a semi-classical approach. These are less reliable, 
especially well below threshold, but have an intuitive 'classical' interpretation in terms of the incoming vacuum 
fluctuations. 



A. Wigner Representation 



In the semi-classical theory, the hierarchy of the stochastic equations given earher can be written, in the frequency 
domain, as follows: 
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FIG. 8: Nonlinear squeezing spectrum with g'' = 0.001, 7r = 0.5 and [i = 0.9. The dashed line represents the analytical result 
and the noisy line the stochastic simulation. 




FIG. 9: Nonlinear squeezing spectrum with = 0.001, 7r = 0.01 and jj, = 0.93. The dashed line represents the analytical 
result and the noisy line the stochastic simulation. This is the driving field for optimum entanglement at zero frequency. 



• First order 



41) (fi) = pl^^ 



-ifl + 7r) 



-in + 1 - 



-iCt + 1 + m) 



-iVt + 1 - 



• Second order 
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• Third order (signal and idler fields) 
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B. Squeezing Correlation spectrum 



The spectrum of the squeezed quadrature, for instance, is given by 

(y(ni)y+(n2)) = g^{y^^\n^)y+^^\n2))+g''{{y^^\nr)y+^-'\n2)) 
+ (y(i)(ni)y+(^Hn2)) + (y(3)(r!i)y+(i)(n2))} + • • ■ 
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The lowest order contribution turns out to be 

Similarly, for the amplified quadrature, to the lowest order, we have 
For the pump quadratures, there is no squeezing, to the lowest order: 



The next contribution to the squeezed quadrature are 

1 - M + 7r 



1 + ^ + 7r 



+ 



(1 + m) [nl + ii + fi + jr)'] 

and 

[^I + (1 + /^)^] 
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+ 



(1+M) [f^2 + (l + A* + 7.)2] 



(l + A. + ^,)»2^(l,^2)(i_^^^^) ^ 

(l-M)[1^2 + (l-M + 7r)2] J ■ ^ ' ' 

This, in turn, gives the following expression for the external squeezing spectrum, obtained by including both internal 
fields and the correlated refiected vacuum noise : 

^ > f22 + (1 + + + (1 + ^)2]2 \ ^^(1 _ ^2) 

_^ [(1 -M)(1 -M + 7r) - + (1 - M + 7r) (1 + + + M^) 



(l-A^) [fi2 + (l-/z + 7.)2] 

[(1 + /Z)(l + /X + 7.) + 2^x2] 02 ^ ^ ^ ^ (1 + 3^ + ^2 _ ^3-) 

(1 + m) + (1 + /i + 7021 ■ ^ ^ 



It is interesting to note that this spectrum is quite different from that given by positive P-rcprcsentation when 
jj, ^ 0. However near the threshold, that is in the limit /x — > 1, the two results show close agreement. This means 
that even when the pump is off, the semi-classical theory gives a distorted vacuum spectrum due to the presence of 
the nonhnear crystal. This happens because in this theory the vacuum fluctuations are taken as real, and then two 
vacuum modes can interact inside the crystal as real fields. In the limit of 7,- 0, the two spectra again become 
compatible, as the semi-classical theory decouples the fundamental mode from its vacuum input in this limit. In the 
case of threshold fiuctuations, we can interpret the agreement as due to the fact that in this region large numbers of 
photon numbers involved - which means that the truncation approximation used in the semi-classical approximation 
is fairly reliable. 
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C. Triple spectral correlation 

For the triple spectral correlation function 

the term proportional to vanishes, and the result, to the lowest non-trivial order is found to be 



(7.12) 



2tt 



1 



i-iiis + ^rmi + (1 - fi^mi + (1 + 

7r 



(7.13) 



As before, the essential difference between quantum and semi-classical theories is that the former gives a zero 
spectrum in the absence of a driving field while the latter, due to the real character of the vacuum field, gives a non 
zero correlation. 

VIII. CRITICAL PERTURBATION THEORY 

As we have seen, the perturbative corrections diverge at the critical point (/i — 1) and a new approach is called for 
to investigate the neighborhood of the threshold. To this end we define new scaled quadratures variables, and use a 
different expansion [23| valid around the critical region. The new pump mode variable xq now corresponds to the real 
scaled depletion in the pump mode amplitude, relative to the undepleted value at the critical point. The signal-idler 
quadrature variables x , a;"*" now describe the critical fiuctuations scaled to be of order 1 at the threshold. 



A. Positive-P Representation 



We scale the quadratures as 



Xo 



7 



x+ = 



, yo ^ 

y+ = Y + 



Yo 



and define also a new scaled time and driving field 



9 V^c 



V = 
T = jgt 

In terms of these variables, the equations in positive-P become 

gdxQ = --ir [xn -2ri + xx'^ - gyy^] dr 
gdya = -7r [va + xy+ + yx^] dr 

dx = - {xqx + gyoy) dr + dw^iir) 



i.2) 



gdy 



+ 2 (xyo - yxo) dr + dwyiir) 
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= - {xqx^ + gyay^) dr + dwx2{T) 



gdy^ 



-2y+ + 2 (^^J/o - y^a^o) dT + dwy2{T) 



i.3) 



The Gaussian white noise sources in these equations are no longer uncorrelated and have the follow the properties 

{dwxidwx2) = 2(1 + |xo^ rfr 



{dwyidwy2) = -2g (^1 + '-xoj dr 
{dwxidwy2) = {dwx2dwyi) = g^yodr 



(8.4) 



We now develop a perturbation theory valid near the threshold. The first set of equations is obtained by neglecting 
all terms of order g or greater on the right sides of the two sets of equations given above 



gdx^^^ 
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(0) 



-(0)^(0) 



dr + dw^2 



V2 



A significant feature of these equations is that the quadratures y^^\ and y^^°\ can be worked out without reference 
to any of other variables, and they give zero noise in the external quadrature at zero frequency. Coupling between 
variables appears in high orders expansion and generates the critical fluctuations in the squeezed quadrature. 

We now consider what happens at or near the classical threshold 77 = 0. In a model where the sub-harmonic 
generation does not cause the pump mode to deplete, we would have x'^^ — 277, and at threshold the critical fluctuations 
in X and x~^ would diffuse outward without any bound. When depletion is included, the critical fluctuations in these 
quadratures are flnite, but very slowly varying compared to those in the other variables. The pump fleld can therefore 
be adiabatically ehminated to first order in the expansion. 

Near threshold {grj <C 1) the decay term in the un-squeezed quadrature x and a;+ is roughly —xq, which is of order 
1. The pump mode will be depleted, so xq must be negative in order for this to be stable. The scaled pump field decay 
is Jr/g, and the squeezed quadrature decay is of order 1/g. If 7^ is much larger than g, it is possible to adiabatically 
ehminate both the pump amplitude and the squeezed quadrature in the equations for the large critical fiuctuations 
X and x^. Since we are taking the limit of small g, we shall assume that this is possible to zero-th order in the 
asymptotic expansion. In the adiabatic elimination, we must solve for the steady state values of the pump xq, given 
an instantaneous first order critical fiuctuation x and x~^ . To leading (zeroth) order this gives 



2ry- 



(8.6) 



Substituting in the equations for x and a;+ , we find that 



After the following change of variables 



r]x 
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r]x 



+ (0) 



+ (0) 



dr + dw^i 



.(0) 



dr + dw^^ 
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the equations 118. 7|l can be put in the form 



X = -ryx- -x(x • x) +Cit) 



(8.9) 



28 



where x is a two component vector whose elements are x+ and X- . 

It is possible to write the Fokker-Planck equation for the probability density t), and look for the equilib- 

rium distribution is of the form -P(x) = Nexp[—U{-x)], where U{-x) is a potential function given by 



1 2 
C/(x) = yyx • X + - (x • x) 

The variance of the critical fluctuations at the critical point, 77 = 0, is given by 

- 2^iil^ = 1.1284 

^ ' r(i/2) 



(8.10) 



(8.11) 



B. Critical Squeezing in positive-P Representation 

We can now flnd the steady state variance of the squeezed quadrature at threshold. Because the fluctuations in the 
squeezed quadrature are very small, we must work to higher order in the asymptotic expansion to obtain a non trivial 
result. To achieve this, it is most useful to introduce equations in the higher order moments yy^ and z = x^y + xy^ . 
The corresponding stochastic equations are derived using Ito rules for the variable changes, so that 



1 + 2yy+ + I [ 2:0 + xoyy+ - ^y^z 



dr + ydwy2 + y^dwyi 



gd{yy+) = -2 
gdz = 

Taking the expectation value at the steady-state {d{yy^)) = 0, we get the flrst order correction 



— 2z + '-ya (2xx'^ + 2gyy^ + Ag) dr + xdwy2 + x^dwyi + gydwx2 + gy^dwxi 



{yy+)^'^ = -^^{{i + yy^)xo~^yoz)^'^ 



The first term in the above expression gives the result 

1 



((l + yy+) :eo)(°) 

For the second term we must write the correlation from the following equation 

gd{yoz) = - [(2 + jr)yoz + lrZ^\ dr + 0{g) + noise 

and then we get 

{yoz)^'^ = -ij^{z'r = -j^{{^^y + ^y+)y"^ = ^ 
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So, finally we get, to first order, 
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(8.12) 
(8.13) 

(8.14) 
(8.15) 
(8.16) 



(8.17) 



This result shows that the best squeezing, in the overall moment, for the intra-cavity combined mode quadrature 
occurs just above threshold, in much the same way as in the degenerate OPO [2ll |. 



C. Wigner Representation 



As in the positive-P equations, we define new scaled quadratures variables to avoid divergences at the critical point 



Xo 



, yo = 

y = Y 

y' 



5.18) 
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In these new variables, the stochastic equations in the Wigner representation are ) 



gdxo 
gdyo 

dx 
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'^^^ = 7i i^ox^ + s/gyoy^) dr + dwx2{T 
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dr + dwyi{T) 



dr + dwy2 (t) 



(8.19) 



Here we use the same notation for scaled time and driving field as in the positive-P case. The noise correlation are 
given by 



{dwxodwyo) = A'y^gdr 
{dwxidwx2) = 2dT 
{dwyidwy2) = 2gdT 



(8.20) 



To develop a perturbation scheme, we define the zero order approximation to be the one in which terms of order 
and greater than are neglected in the set of equations above 
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(8.21) 



It is worth noting that this set of equations, though having the same structure as that in the positive-P case, has 
differences in the correlations of the noise terms. On adiabatic ehmination of the pump and substituting this result 
into x^ and x~^^ we find the same equations as in the positive-P representation, since to zero-th order the correlation 
noise in both theories is identical. 



D. Critical squeezing in Wigner representation 

Now we proceed to calculate {yy~^) at threshold using the Wigner representation. Using the It6 rules we get 

gd{yy+) = 2- Ayy+ + ^VoZ - ^2yy+xo + dwyi + dwy2 (8.22) 

where we have defined z = yx~^ + y~^x, which obey the following equation 

gdz = -2z + ,/gyoxx~^ + g^/gyoyy^ + x'^dwyi + gydwx2 + xdwy2 + gy^dwxi (8.23) 

The squeezing variance at threshold in the steady state is obtained from the above equation taking expectation 
values 

= \ + ^{yoz) - l{xoyy+) (8.24) 

The last term of the above equation can be written as 

f(4°^)(2/y+)^°^ = f (8.25) 
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and the equation (55) gives the result 

{yoz^ = -Vff^^^ + Vg ^"°^ ^ (8.26) 

Using the results derived from the zero order equations 

^ 2(xWa;+(°))(y(")2;+(")) (8.27) 

we finally get 

This result is exactly the same as obtained in positive P-representation. We can infer that dropping third order 
terms in the Wigner phase space equation does not have any direct consequence for the near threshold analysis of 
entanglement to this order of approximation. This is to be contrasted with the situation far below threshold, where 
there are large differences in the nonlinear contributions, indicating a failure of the truncated (hidden- variable) Wigner 
theory. 

The change in behavior has a simple mathematical origin. Far below threshold, the signal/idler photon numbers are 
small, which leads to a failure of the truncation approximation when using the semi-classical method. At the critical 
point, photon numbers in all modes are relatively large, so the truncation approximation has less severe consequences. 



IX. CONCLUSIONS 



We have calculated the effects of nonlinear quantum fiuctuations in a nondegenerate parametric oscillator, both 
below and at the classical threshold, using stochastic equations that follow from the positive P-representation. The 
analytical results thus obtained are compared with exact numerical simulations. The spectral entanglement and 
squeezing in the output fields is maximized just below threshold. This may be useful, for example, in cryptographic 
applications!^. We find that at the critical point (p, = 1), the scaling behavior is quite different to the behavior 
below threshold, and must be calculated by using an asymptotic perturbation theory, valid at the threshold itself. The 
total intra-cavity squeezing and entanglement moment is actually minimized at a driving field just above threshold. 
This behavior was confirmed in our simulations. This apparent paradox can be attributed to the fact that the critical 
fluctuations mostly tend to broaden the squeezing spectrum, which has a strong effect at zero-frequency but does not 
diminish the total squeezing moment, integrated over all frequencies. 

A similar analysis was carried out within the framework of the semi-classical theory arising from a truncation to a 
Fokker-Planck form of the evolution equation in the Wigner representation. Here, we found that well below threshold, 
while the linear terms agreed with full quantum calculation, the nonlinear corrections tend to disagree, especially for 
low sub-harmonic losses. However, at the critical point, the situation changes. Here, where the dominant terms are 
nonlinear, we find excellent agreement between the two methods. While quantum fiuctuations are indeed large at 
the critical point, it appears that an equally acceptable interpretation of the observed noise characteristics near the 
critical point exists via a semi-classical model, which is essentially a kind of hidden- variable theory. 

Our main result is that entanglement, EPR correlations and squeezing are optimized very near threshold. At the 
same time, the semi-classical Wigner approximation can give an excellent description of the squeezing and entangle- 
ment fiuctuations near threshold. On the other hand, some highly nonclassical signatures of quantum correlations 
occur in the higher-order correlations, which are not described by the semi-classical approach. Surprisingly, these 
nonclassical and non-Gaussian signatures only occur well below threshold, where one might have expected the usual 
linearized analysis to be applicable. 

This suggests that experimental tests of the present theory may be carried out either near threshold - where the 
largest effects will be observed in the enhanced critical fluctuations of the unsqueezed quadrature - or well below 
threshold, where nonclassical triple correlations are predicted. 
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